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Abstract 



We consider a three-dimensional lattice model consisting of layers of ver- 
tex models coupled with interlayer interactions. For a particular non-trivial 
interlayer interaction between charge-conserving vertex models and using a 
transfer matrix approach, we show that the eigenvalues and eigenvectors of 
the transfer matrix are related to those of the two-dimensional vertex model. 
The result is applied to analyze the phase transitions in a realistic three- 
dimensional dimer system. 
05.50.+q 
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An outstanding unsolved problem in the statistical mechanics of critical phenomena is 
the pursuit of exact solutions for realistic three-dimensional (3-D) systems. While a large 
number of 2D systems have yielded to analyses only a limited number of 3D systems 
have been solved. They include the 3D Ising model solved by Suzuki 0, the Zamolodchikov 
model [[J solved by Baxter || and its more recent iV-state extension by Bazhanov and 
Baxter ||. However, these models invariably suffer defects in one way or another: The 
Suzuki model turns out to be a 2D system in disguise, while the Zamolodchikov model and 
its extension involve unphysical negative Boltzmann weights. Similarly, a continuous string 
model in general d dimensions solved by two of us || also involves negative weights. The 
solution of realistic physical 3D models has remained very much an open problem. 

One approach toward solving realistic 3D models is to build from 2D systems. Indeed, 
such an approach has been suggested 0] and applied with some success recently . How- 
ever, in these considerations much attention has been placed to the algebraic structure of 
the transfer matrix and the associated Yang-Baxter equation, to the extent of masking the 
elegance of the solution. Here, we present a more general formulation, and deduce its solu- 
tion via an alternate and yet much simpler consideration. The result is applied to analyze 
a realistic layered 3D dimer system. 

Consider a simple-cubic lattice £ of size K x M x N with periodic boundary conditions. 
Regard £ as consisting of K copies of square lattices, of M rows and N columns each and 
stacked together as shown in Fig. 1. For simplicity, we shall speak of the square lattices 
as "layers" of £. Label sites of £ by indices {m,j, k}, with 1 < m < M, 1 < j < N and 
1 < k < K. Within each layer of £ define a 2D g-state vertex model whose lattice edges 
can be in q distinct states. Label the state of the horizonal (vertical) edge incident at the 
site {m,j,k} in the direction of, say, decreasing {m,j} by a m jk (Pmjk)- It is convenient at 
times to suppress the subscripts m and/or k by adopting the notation 

and similarly for the a's. Associate vertex weight W m jk to site {m,j, k} which is a function 
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of the configuration C md>k = {a mjk , (3 mjk , a m>j+ i, k , /3 m+ i d>k } -> {a j} fy, a j+1 , /3-} of the four 
edges incident at the site {m,j, k}. Let the {m, j} sites of two adjacent layers k and k + 1 
interact with a Boltzmann factor B m j k which in the most general case is a function of the 
configurations C m j tk and C m j,fc+i of the eight edges incident to the two sites. Then, the 
problem at hand is the evaluation of the partition function 

K m n , v 
Zmnk — E E TT IT TT( B m j k W m j k J (2) 

Ctmjk Pmjk k = i m=i j = i ^ 

where the summations are taken over all edge states a m j k and /3 mj -fc, and the per-site "free 
energy" for any K 

f K = K' 1 JunjMN)- 1 In Z MNK . (3) 

The transfer matrix: The partition function (Q) can be evaluated by applying a transfer 
matrix in the vertical direction. In a horizontal cross section of £ there are NK vertical 
edges. Let {f3 m } = {fl m j k \l < j < N, 1 < < K}, 1 < m < M denote the states of these 
NK vertical edges, and define a 2 NK x 2 NK matrix T with elements 

K N , s 

T({/? m }, {/W) = E II II f^miJfcWmifc), 1 < m < M. (4) 

arnjfc fc = li = l ^ ^ 

Then one has 

AI 

Zmnk = E IT ^({/^m}, {/?m+i}) 
= Tr 

- A^fax (5) 

where A max is the largest eigenvalue of T. 

It is clear that we need to restrict considerations to models which are soluble when 
the interlayer interaction is absent, or B m j k = 1. This leads us to build 3D systems from 
soluble 2D models. It is also clear that the interlayer interaction B m j k should be such that 
the overall interlayer factor Y[ m ,j, k B m j k can be conveniently treated. For this purpose we 
restrict considerations to 2D charge- conserving models. 



For definiteness let the labels a m jk and /3 m jk take on a set I of q integral values. For 
example, one can take I = {+1, —1} for q = 2 and / = {+1, 0, —1} for q = 3. A 2D vertex 
model is charge- conserving if its vertex weights are non-vanishing only when 

a j + Pj = (Xj+x + P'j (charge conservation) (6) 



holds at all sites. Examples of charge conserving models are the q = 2 ice- rule models ||10|| . 
the g-state string model |[L1|| , the q = 3 Izergin-Korepin model ]nj, and others |13 



A direct consequence of the charge-conserving rule (^) is deduced by summing (|B|) from 
j = 1 to j — N. This yields 

i JV I N 

3=1 J=l 

showing that the quantity — 1 < yfc < 1 is independent of m. (For ice-rule models this fact 
is well-known.) 

Next one introduces the interlayer interaction [0-0 

B mjk = exp(h(ajj3j - oij+ifi'M. (8) 

Since the negation of h corresponds to a reversal of the layer numberings, without loss of 
generality we can take h > 0. We now show quite generally that the interlayer interaction 
(H) leads to a considerable simplification of the transfer matrix. Consider first the product 

N , N s 

II B mjk = exp ( hJ2( a jPj - Oij+iP'j) ) • (9) 

3=1 V 3=1 J 

Summing over (^), or ctj + fa = a i+ i + (3'^ for i = {1, j — 1} and i = { j + 1, iV} for the layer 
k + 1, one obtains, respectively, the identities 

j-i 

=«i + E(A-A0, J = 2,3---,iV 
1=1 

TV 

fii+i = fii- E(A-A'), i = i,2,-..,iv-i (10) 

i=3+l 

where we have used ajv+i = «i- Substituting (|T0|) into (§) and making use of the identity 
J2jL2 Si=i = S^i 1 SiLj+i m the first summation in (|j), one arrives after a little algebra at 



n B mjk = ex P (Nh( ai y - a lV ) + NhtftfJ) - f(P'J'))), (11) 

where /(/3,/3) = J2f=2 Y$Zi PiPj- The numerical factor /(/3, /3), which is defined for each 
fixed m, is cancelled in the further product 

M N M , v 

II II ^rryfc = II ex P ( Nh(a m ^ k y k+1 - a mi i jfc+ i?/ fc )j- (12) 

m=l j'=l m=l 

As a result, only the conserved quantities y k and the state a m ,i,k °f the extremites of a row of 
horizontal edges appear in the product (|TJ)- This leads us to rewrite the partition function 
(1) as 

Z MNK = Tr (T cfr ) M (13) 
where T eff is a matrix with elements 

K f N \ 

T eff ({/U, {/W) = E II (e^-.^fe+i-^-i) 7J Wmjfc ). (14) 

"mjt fc=l ^ j = l ' 

The problem is now reduced to one of finding the largest eigenvalue of T eff . In fact, expression 
( p|) shows that T efT is related to T by a similarity transformation T eff = S T S _1 where S 
is diagonal. It follows that T and T eff have the same eigenvalues, and their eigenvectors are 
related. The task is now considerably simpler since one needs only to keep track of the 2D 
system. The problem is solved if the eigenvalues of the transfer matrix for the 2D vertex 
model can be evaluated for fixed y k and a m ^ jk . 

The ice-rule model: To illustrate the usefulness of this formulation, we now apply it to 
layers of ice-rule model with vertex weights {ui,U2, (for standard notations relevant 

to present discussions see, for example, |10|). Let a = +1 (—1) denote arrows pointing 
toward right (left), and f3 — +1 (—1) arrows pointing up (down). Then one verifies that the 
charge- conserving condition (§) is satisfied with y k = 1 — 2n k /N, where n k is the number 
of down arrows in a row of vertical edges in the kth layer. Introducing next the interlayer 
interaction (§), the eigenvalues of the matrix (14) are obtained by applying a global Bethe 
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ansatz consisting of the usual Bethe ansatz for each layer. The algebra is straightforward 
and one obtains 

K 

Zmnk ~ max TT [A R (n k ) + A L (n fc )] M , (15) 

Kni.<N 
- K ~ k=l 



with 



A R (n k ) = e -2M»*+i-»*-i) w f -»* Yl 



\ (k) 
j=l v u 4 - UJiZ) 



3=1 



(k) 



where A# (Aj,) refers to the eigenvalue for a m ^ k = +1 (—1) |L0| and, for each 1 < k < K, 
the n k complex numbers Zj k \ j = 1, 2, . . . , n k are the solutions of the Bethe ansatz equations 

e ±h{n k+l -n u _ x )/ z {k)^N = (_-\\n k +lj\( B{ z ii z j 
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ir +1 n(^() ; J = ^,...,n k (17) 



1 = 1 "I) 

where B(z, z') = co> 2 ^4 + ujiuj^zz' — (uiU2 + uj 3 u 4 — uj§lo§)z' . Note that the Bethe ansatz 
equation ([T7|), which is obtained by imposing the cancellation of unwanted terms in the 
Bethe ansatz solution, differs from its usual form (see, for example, WW) in the inclusion 



of the exponential factor involving h. Various special forms of this solution has been given 
previously |7| 0- 

A dimer system with interlayer interactions: We now consider a 3D lattice model con- 
sisting of layers of honeycomb dimer lattices. The dimers, which carry weights u, v, w along 
the three honeycomb edge directions, are close packed within each layer and, in addition, 
interact between layers. For two dimers incident at the same {m, j} site in adjacent layers, 
the interaction energy is given in Table I. 

The 2D honeycomb dimer system can be formulated as a five-vertex model, namely, an 
ice-rule model with the weights [ |i4] , [r"5[ 



{d>i, LU 2 , ^3, ^4, k> 5 , LU e } = {0, W, V, U, y/UV, ^/iw} . (18) 

The 5-vertex model is defined on a square lattice of size M x N mapping to an honeycomb 



lattice of 2MN sites | 14)J15 |. The mapping is such that the edge state a = +1 ((3 = +1) 



corresponds to the presence of a v (u) dimer. It can then be verified that the interlayer 
interaction given in Table I can be written precisely in the form of (^j) ||16| , and therefore 
we can use the ice-rule model results. 

Substituting (|IBD into (0), one obtains 

K nk /111 V „A M 

Zmnk - u max [] II (~ + \ f ) ^ 

l<n k <N ^^Vu U J J 

with the Bethe ansatz solution 

z f = e «>i e Mvk + i-Vk-i) i j = i ; 2 , ■ ■ ■ , n k (20) 

where e t6j are n k distinct A^th roots of (— l) nfc+1 . For a given n k , the factor inside the 
parentheses in ( |19"D attends its maximum if the 9/s lie on an arc crossing the positive real 
axis and extending from — 7r(l — yk)/2 to 7r(l — t/fc)/2. Using (|3]) this leads to the per-site 
free energy 

f K = Inn + max - V — f^' 2 + ^ e 2My ft+ i-%-x) e i^^. (2 i) 

This is our main result. 

We have carried out analytic as well as numerical analyses of the free energy ( |2"T| ) for 
K = 3 x integer. Here we summarize the findings. For h = 0, the layers are decoupled and 
the property of the system is the same as that of the 2D system |L4|,[ll|. For large h, it is 
readily seen from Table I that the energetically preferred state is one in which each layer is 
occupied by one kind of dimers, u, v, or w, and the layers are ordered in the sequence of 
{w, v, it, w, v , u, ■ ■ ■}. It is also clear that for large u, v, or w, the system is also frozen with 
complete ordering of u, v, or w dimers. These orderings are referred to as the H, U, V, and 
W phases, respectively. The system can also be in two other phases. A Y phase in which 
all layers have the same value of = y determined straightforwardly by maximizing (pT|). 
or 

w 2 + v 2 + 2wv cosf|(l-?/)J =m 2 , (22) 
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and an I phase which is the H phase with any of the w, v, or u layers replaced by layers with 
yk = y- If the v layer is replaced by a y layer so that the ordering is {w, y, u, w, y,u, ■ ■ ■}, 
for example, then y is given by ( p2|) with v replaced by ve Ah . 

The phase diagram is found to be symmetric in w, v, and u. It is then convenient to 
plot the phase diagram using the coordinates 

X = \n(v/w) Y = (v^r 1 \n(vw/u 2 ) (23) 

so that any interchange of the three variables w, v, and u corresponds to a 120° rotation 
in the {X, Y } plane. The phase diagram for h < h Q = 0.2422995... is the same as in Fig. 
2a but without the H regime. Increasing the value of h one finds the H phase appear in 
h < h < hi = 0.2552479... as shown in Fig. 2a. At h = hi the / phase appears (Fig. 
2b), with its region extending to infinity when h reaches h 2 = (ln3)/4 = 0.2746531... (Fig. 
2c). When h reaches h% = 0.3816955... and higher, the Y phase disappears completely as 
shown in Fig. 2d. All transitions are found to be of first-order except the transitions between 
the {U, V, W} and Y phases, and between the I and H phases, which are found to be of 
second-order with a square-root divergence in the specific heat. 

In summary, we have presented the formulation of a general 3D lattice model, and applied 
it to solve a realistic dimer system. The analysis can be extended to include dimer-dimer 
interactions within each layer ||14|| , and further-neighbor interlayer interactions. Details of 



the present and further analyses will be presented elsewhere. Work has been supported in 
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FIGURES 

FIG. 1. A three-dimensional lattice model consisting of layered vertex models. 

FIG. 2. Phase diagrams of the dimer system, (a) ho < h < hi, (b) hi < h < hi, (c) 
h 2 <h< h 3 , (d) h > h 3 . 
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TABLES 

TABLE I. Interaction energy between two dimers incident at the same {m, j} site of adjacent 
layers. For example, a v dimer in the kth layer interacts with a u dimer in the (k + l)th layer with 
an energy eh. Here e = +1 (—1) if the site is in sublattice A (B). 

layer k — ► k + 1 w v u 

w —h h 

v h eh 

u —h eh 
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